home *** CD-ROM | disk | FTP | other *** search
/ Aminet 52 / Aminet 52 (2002)(GTI - Schatztruhe)[!][Dec 2002].iso / Aminet / misc / emu / Apex-src.lha / TRANSCND.XPL < prev    next >
Text File  |  2001-09-30  |  9KB  |  425 lines
  1. \TRANSCND.XPL    OCT-28-86
  2.  
  3. code
  4. ABS=0,        RAN=1,        REM=2,        RESERVE=3,
  5. SWAP=4,        EXTEND=5,    RESTART=6,    CHIN=7,
  6. CHOUT=8,    CRLF=9,        INTIN=10,    INTOUT=11,
  7. TEXT=12,    OPENI=13,    OPENO=14,    CLOSE=15,
  8. FIX=50,        HEXOUT=27;
  9. code real    RLRES=46,    FLOAT=49,    RLABS=51;
  10.  
  11. int    II, JJ, M, N, MDIGIT;
  12. real    XX, YY, AA, BB;
  13.  
  14.  
  15. \def    PI= 3.1415926535897932,        \\PI
  16. \    PI2= 6.2831853071795864,    \\PI *2
  17. \    HALFPI= 1.5707963267948966,    \\PI /2
  18. \    HALFPI3= 4.7123889803846898;    \\PI /2 *3
  19.  
  20. real    PI, PI2, HALFPI, HALFPI3;
  21.  
  22.  
  23.  
  24. func real SQRT(X);    \SQUARE ROOT FUNCTION
  25. real    X;
  26. real    G1, G2;
  27. int    EXP, I;
  28. addr    A, A2;
  29. begin
  30. if X < 0.0 then X:= -X;        \ERROR: SQUARE ROOT OF A NEGATIVE NUMBER
  31. A:= addr X;
  32. EXP:= SWAP(A(0)) + A(1);
  33. EXP:= EXP /2 + $1FF8;
  34. A2:= addr G1;
  35. G1:= 2.0;
  36. A2(0):= SWAP(EXP); A2(1):= EXP;
  37. for I:= 0,6 do
  38.     begin
  39.     G2:= X /G1;
  40.     G1:= (G1 + G2);
  41.     G1:= G1 /2.0;
  42.     if G1 = G2 then I:= 1000;
  43.     end;
  44. return G1;
  45. end;    \SQRT
  46.  
  47.  
  48.  
  49. func real MOD(X, Y);    \MODULO FUNCTION
  50. \E.G: MOD(10.0, 3.0) = 1.0; MOD(-12.3, 1.0) = -0.3
  51. \WARNING: THIS ROUTINE LOSES PRECISION AS RESULT APPROACHES 0.0
  52. \ ALSO, ATTEMPTING TO GET MODULOS BEYOND 32 BITS OF PRECISION CAUSES
  53. \ A FIX OVERFLOW.
  54. real    X, Y;
  55. real    Z;
  56. begin
  57. Z:= X /Y -0.5;
  58. return X - FLOAT(FIX(Z)) *Y;
  59. end;    \MOD
  60.  
  61.  
  62.  
  63. func real POLLY(X, P, N);
  64. \EVALUATE X USING A POLYNOMIAL EXPRESSION OF THE FORM:
  65. \ = P0 + P1 *X^2 + P2 *X^4 + P3 *X^6
  66. \ = (((P3 *X^2) + P2) *X^2 + P1) *X^2 + P0
  67. real    X,    \NUMBER TO APPROXIMATE
  68.     P;    \TABLE OF COEFFICIENTS
  69. int    N;    \NUMBER OF TERMS
  70. int    I;    \LOOP COUNTER
  71. real    X2;    \X^2, X SQUARED
  72. begin
  73. X2:= X *X;
  74. X:= P(N-1);
  75. for I:= -(N-2), 0 do
  76.     X:= X *X2 + P(-I);
  77. return X;
  78. end;    \POLLY
  79.  
  80.  
  81.  
  82. func real COS(X);    \RETURN THE COSINE OF X
  83. real    X;
  84. real    P;
  85. begin
  86. \REDUCE RANGE TO 0 <= X < HALFPI
  87. X:= RLABS(X);
  88. if X>PI2 then X:= MOD(X,PI2);
  89.  
  90. P:=[     0.9999999999999999960897E+00,
  91.     -0.49999999999999974308584E+00,
  92.      0.4166666666666387895916E-01,
  93.     -0.138888888887731721151E-02,
  94.      0.24801587277443938629E-04,
  95.     -0.275573163935346178E-06,
  96.      0.20876561960112253E-08,
  97.     -0.114629048993344E-10,
  98.      0.46090073769E-13];
  99.  
  100. if X<HALFPI then                \QUADRANT 1
  101.     return if X<1E-10 then 1.0 else POLLY(X, P, 9);
  102. if X<HALFPI3 then return -POLLY(X-PI, P, 9);    \QUADRANTS 2 OR 3
  103. return POLLY(PI2-X, P, 9);            \QUADRANT 4
  104. end;    \COS
  105.  
  106.  
  107.  
  108. func real SIN(X);    \RETURN THE SINE OF X
  109. real    X;
  110. begin
  111. return if RLABS(X) < 1E-5 then X else COS(HALFPI -X);
  112. end;    \SIN
  113.  
  114.  
  115.  
  116. func real TAN(X);    \RETURN THE TANGENT OF X
  117. real    X;
  118. real    Y, Z;
  119. begin
  120. Y:= SIN(X);
  121. Z:= COS(X);
  122. return Y /Z;
  123. end;    \TAN
  124.  
  125.  
  126.  
  127. func real ATAN(X);    \ARC-TANGENT FUNCTION
  128. real    X;
  129. real    P,
  130.     ATAN25,
  131.     ATAN75;
  132.  
  133.  
  134.     func real ATANY(X);    \ARC-TANGENT FUNCTION FOR 0.0 >= X < 1.0
  135.     real    X;
  136.     real    Z;
  137.     begin
  138.     if X >= 0.5 then
  139.         begin
  140.         Z:= (X -0.75) /(1.0 + X *0.75);
  141.         X:= POLLY(Z, P, 9) *Z + ATAN75;
  142.         end
  143.     else    begin
  144.         Z:= (X -0.25) /(1.0 + X *0.25);
  145.         X:= POLLY(Z, P, 9) *Z + ATAN25;
  146.         end;
  147.     return X;
  148.     end;    \ATANY
  149.  
  150.  
  151.     func real ATANX(X);    \ARC-TANGENT FUNCTION FOR POSITIVE X
  152.     real    X;
  153.     begin
  154.     return if X >= 1.0 then HALFPI -ATANY(1.0/X) else ATANY(X);
  155.     end;    \ATANX
  156.  
  157.  
  158. begin
  159. ATAN25:= 0.244978663126864154;
  160. ATAN75:= 0.643501108793284386;
  161. P:= [     0.9999999999999999849899E+00,
  162.     -0.333333333333299308717E+00,
  163.      0.1999999999872944792E+00,
  164.     -0.142857141028255452E+00,
  165.      0.11111097898051048E+00,
  166.     -0.909037114191074E-01,
  167.      0.767936869066E-01,
  168.     -0.6483193510303E-01,
  169.      0.443895157187E-01];
  170.  
  171. return if X < 0.0 then -ATANX(-X) else ATANX(X);
  172. end;    \ATAN
  173.  
  174.  
  175.  
  176. func real ATAN2(Y,X);    \RETURN THE ARCTANGENT OF Y/X
  177. real    Y, X;
  178. begin
  179. if X = 0.0 then return HALFPI *Y /RLABS(Y);
  180. if X > 0.0 then
  181.     return ATAN(Y/X)
  182. else    return if Y >= 0.0 then ATAN(Y/X) + PI else ATAN(Y/X) - PI;
  183. end;    \ATAN2
  184.  
  185.  
  186.  
  187. func real ASIN(X);    \ARC-SINE FUNCTION
  188. \WARNING: INACCURATE AS RLABS(X) APPROACHES 1.0.
  189. real    X;
  190. real    Z;
  191. begin
  192. Z:= SQRT(1.0 - X *X);
  193. return ATAN(X /Z);
  194. end;    \ASIN
  195.  
  196.  
  197.  
  198. func real ACOS(X);    \ARC-COSINE FUNCTION
  199. \WARNING: INACCURATE AS RLABS(X) APPROACHES 1.0.
  200. real    X;
  201. real    Z;
  202. begin
  203. return -ASIN(X) + HALFPI;
  204. end;    \ACOS
  205.  
  206.  
  207.  
  208. proc    FORMAT(M1, N1);    \SET FORMAT PARAMETERS FOR RLOUT
  209. int    M1, N1;
  210. begin
  211. MDIGIT:= M1;
  212. N:= N1;
  213. end;    \FORMAT
  214.  
  215.  
  216.  
  217. proc    RLOUT(DEV, X);
  218. \Output the real number X to the specified device.
  219. \Other inputs: M, N.
  220. int    DEV;
  221. real    X;
  222. real    SX, RND, HALF, ONE, TEN;
  223. int    I, K, L, NEG;
  224. def    SIGFIGS =15;    \Maximum number of decimal digits
  225.  
  226.  
  227.     proc    DIGITOUT;
  228.     int    DIGIT;
  229.     begin
  230.     for I:= 1, K do
  231.         begin
  232.         if L > 0 then
  233.             begin
  234.             X:= X *TEN;
  235.             DIGIT:= FIX(X -HALF);
  236.             CHOUT(DEV, DIGIT +^0);
  237.             X:= X -FLOAT(DIGIT);
  238.             L:= L -1;
  239.             end
  240.         else    CHOUT(DEV,^0);
  241.         end;
  242.     end;    \DIGITOUT
  243.  
  244.  
  245. begin
  246. TEN:= FLOAT(10);
  247. ONE:= FLOAT(1);
  248. HALF:= ONE /FLOAT(2);
  249. if X < FLOAT(0) then [X:= -X; NEG:= true] else NEG:= false;
  250.  
  251. K:= 0;
  252. SX:= X;            \SAVE ORIGINAL NUMBER TO DETERMINE LEADING ZERO
  253. if X # FLOAT(0) then
  254.     begin
  255.     while X >= ONE do [X:= X /TEN; K:= K +1];
  256.  
  257.     \ADD IN ROUNDING FACTOR: 0.5 * 10 ^ -(K+N)
  258.     RND:= HALF;
  259.     L:= K +N;
  260.     if L > SIGFIGS then L:= SIGFIGS;
  261.     for I:= 1, L do RND:= RND /TEN;
  262.     X:= X +RND;
  263.  
  264.     if X >= ONE then
  265.         [X:= X /TEN; K:= K +1;    \ADJUST FOR ROUND OVERFLOW
  266.         SX:= TEN];        \FORGET ABOUT LEADING ZERO
  267.     end;
  268.  
  269. \Calculate the number of leading blanks needed:
  270. L:= M -K;
  271. if SX < ONE then L:= L-1;        \LEAVE ROOM FOR LEADING ZERO
  272. for I:= 1, L do CHOUT(DEV,^ );
  273. CHOUT(DEV, if NEG then ^- else ^ );
  274. if SX < ONE then CHOUT(DEV,^0);        \OUTPUT LEADING ZERO, E.G: 0.2
  275.  
  276. L:= SIGFIGS;
  277. DIGITOUT;                \OUTPUT DIGITS IN FRONT OF THE D.P.
  278. if N > 0 then [CHOUT(DEV,^.); K:= N; DIGITOUT];  \OUTPUT DIGITS AFTER D.P.
  279. end;    \RLOUT
  280.  
  281.  
  282.  
  283. proc    RLOUTX(DEV, X);
  284. \Output the real number X to the specified device.
  285. \Other inputs: MDIGIT, N.
  286. int    DEV;
  287. real    X;
  288. int    NEG, EXP;
  289. real    ZERO, ONE, TEN, KILO;
  290.  
  291.  
  292.     proc    EXPOUT;
  293.     begin
  294.     if NEG then X:= -X;
  295.     RLOUT(DEV, X);
  296.     CHOUT(DEV, ^E);
  297.     CHOUT(DEV, if EXP < 0 then ^- else ^+);
  298.     EXP:= ABS(EXP);
  299.     if EXP < 10 then CHOUT(DEV, ^0);
  300.     INTOUT(DEV, EXP);
  301.     end;    \EXPOUT
  302.  
  303.  
  304. begin
  305. if MDIGIT > 1 then [M:= MDIGIT; RLOUT(DEV, X); return];
  306.  
  307. ZERO:= FLOAT(0);
  308. ONE:= FLOAT(1);
  309. TEN:= FLOAT(10);
  310. KILO:= FLOAT(1000);
  311.  
  312. if X < ZERO then [X:= -X; NEG:= true] else NEG:= false;
  313. EXP:= 0;
  314.  
  315. if MDIGIT = 0 then            \SCIENTIFIC NOTATION
  316.     begin                \E.G: 1.2E+23, 1.2E-102, 1.2E+02
  317.     M:= 2;
  318.     if X # ZERO then
  319.         begin
  320.         while X < ONE do [X:= X *TEN; EXP:= EXP -1];
  321.         while X >= TEN do [X:= X /TEN; EXP:= EXP +1];
  322.         end;
  323.     EXPOUT;
  324.     end
  325. else    begin                \ENGINEERING NOTATION
  326.     M:= 4;
  327.     if X # ZERO then
  328.         begin
  329.         while X < ONE do [X:= X *KILO; EXP:= EXP -3];
  330.         while X >= KILO do [X:= X /KILO; EXP:= EXP +3];
  331.         end;
  332.     EXPOUT;
  333.     end;
  334. end;    \RLOUTX
  335.  
  336.  
  337.  
  338. func real RLIN(DEV);
  339. \Read in the ASCII representation of a real number from the specified device
  340. \ and return its value.
  341. int    DEV;    \Input device
  342. int    CH,    \Character
  343.     EX,    \Power-of-ten exponent, total effective value
  344.     N,    \Exponent as specified by input
  345.     NEG,    \Flag: Negative real number
  346.     ENEG,    \Flag: Negative exponent
  347.     DIGIT;    \Flag: Last character is a digit (0 thru 9)
  348. real    X,    \Value of real number
  349.     TEN;    \1.0, Avoids use of real constants which are not as easily
  350.         \ ported from one floating point representation to another.
  351.  
  352.  
  353.     proc    GETCH;        \Get character from input device
  354.     begin
  355.     CH:= CHIN(DEV);
  356.     DIGIT:= CH>=^0 & CH<=^9;    \Is it a digit?
  357.     end;    \GETCH
  358.  
  359.  
  360.     proc    ADDIN;
  361.     begin
  362.     X:= X *TEN + FLOAT(CH -^0);
  363.     end;    \ADDIN
  364.  
  365.  
  366. begin    \RLIN
  367. TEN:= FLOAT(10);
  368. NEG:= false;
  369. loop    begin
  370.     GETCH;            \Ignore any leading garbage
  371.     if CH =^- then NEG:= not NEG;
  372.     if DIGIT then
  373.         begin
  374.         X:= FLOAT(CH -^0);
  375.         loop    begin
  376.             GETCH;
  377.             if not DIGIT then quit;
  378.             ADDIN;
  379.             end;
  380.         quit;
  381.         end;
  382.     if CH=^. then [X:= FLOAT(0); quit];
  383.     end;
  384. EX:= 0;
  385. if CH = ^. then
  386.     loop    begin
  387.         GETCH;
  388.         if not DIGIT then quit;
  389.         ADDIN;
  390.         EX:= EX -1;    \if X gets bigger, the exponent gets smaller
  391.         end;
  392. if CH=^E ! CH=^e then
  393.     begin
  394.     N:=0;
  395.     GETCH;
  396.     if CH = ^- then [ENEG:= true; GETCH] else ENEG:= false;
  397.     if CH = ^+ then GETCH;
  398.     while DIGIT do [N:= N *10 +(CH -^0); GETCH];
  399.     EX:= EX + (if ENEG then -N else N);
  400.     end;
  401. while EX < 0 do [X:= X /TEN; EX:= EX +1];
  402. while EX > 0 do [X:= X *TEN; EX:= EX -1];
  403. return if NEG then -X else X;
  404. end;    \RLIN
  405.  
  406.  
  407.  
  408. begin    \MAIN
  409. \DEFINED CONSTANTS DON'T WORK
  410. PI:= 3.1415926535897932;
  411. PI2:= 6.2831853071795864;
  412. HALFPI:= 1.5707963267948966;
  413. HALFPI3:= 4.7123889803846898;
  414.  
  415. loop    begin
  416.     FORMAT(0,10);
  417.     YY:= RLIN(0);
  418.     XX:= RLIN(0);
  419.     YY:= ATAN2(YY,XX);
  420.     RLOUTX(0,YY); CRLF(0);
  421.     end;
  422. end;    \MAIN
  423. 
  424. loop    begin
  425.     FORMAT(0,10)